Density-based mixing parameter for hybrid functionals 
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A very popular ab-initio scheme to calculate electronic properties in solids is the use of hybrid 
functionals in density functional theory (DFT) that mixes a portion of Fock exchange with DFT 
functionals. In spite of their success, a major problem still remains, related to the use of one 
single mixing parameter for all materials. Guided by physical arguments that connect the mixing 
parameter to the dielectric properties of the solid, and ultimately to its band gap, we propose a 
method to calculate this parameter from the electronic density alone. This method is able to cut 
significantly the error of traditional hybrid functionals for large and small gap materials, while 
retaining a good description of structural properties. Moreover, its implementation is simple and 
leads to a negligible increase of the computational time. 



Density functional theory (DFT) is one of the major 
achievements of theoretical physics in the last decades. It 
is now routinely used to interpret experiments or to pre- 
dict properties of novel materials. The success of DFT 
relies on the Kohn-Sham (KS) scheme and the existence 
of good approximations for the unknown exchange and 
correlation (xc) functional. In the standard KS formula- 
tion the xc potential is local and static. Since the original 
suggestion of the local-density approximation (LDA) [l| , 
a swarm of functionals has been proposed in the liter- 
ature (2). In the ab-initio study of solids, the Perdew, 
Burke and Ernzcrhof [5] (PBE) parametrization of the xc 
functional has been for many years the default choice for 
many applications. A good functional must yield ground 
states properties (like structural parameters), while it is 
expected that the KS gap and true quasiparticle gap dif- 
fer by the derivative discontinuity |3j. Indeed, for semi- 
conductors and insulators PBE yields good structural 
properties and KS band-gap energies that are at best half 
of their experimental value. To obtain both the ground 
state and quasiparticle energies correctly within one and 
the same formalism,one can resort to, e.g., a many-body 
GW calculation [alzl- However, GW is by all measures 
an expensive technique, with a very unfavorable scaling 
with the number of atoms in the unit cell. It is there- 
fore unpractical for the study of band structures of large 
systems and clearly prohibitive regarding total energy 
calculations even for simple realistic systems. 

Much of the computational effort in GW comes from 
the dynamically screened Coulomb interaction W . It has 
therefore been crucial to explore to which extent dynam- 
ical effects are mandatory, or whether non-locality is the 
dominating characteristic. The move from local KS po- 
tential to non-local functionals has first been pushed for- 
ward in Quantum Chemistry, where today the so-called 
hybrid functionals are very popular. These functionals 
mix a fraction a of Fock exchange with a combination of 
LDA and generalized gradient (GGA) functionals. The 
application of hybrid functionals to the solid state had 
a much slower start The situation changed re- 



cently, helped by the wider availability of computer codes 
that support hybrids fiol ] and the continued increase of 
computational power covering the additional cost with 
respect to a local potential. Besides yielding good struc- 
tural properties [12j | hybrids have proved to correct to 
a large extent the band-gap problem [tj [ll|. Another 
landmark came with the introduction of screened hy- 
brids 13]. These functionals lead to faster calculations 
and improved band-gaps, especially for small band-gap 
systems. Furthermore, by screening the Coulomb inter- 
action at large distances, they also give access to metals. 

The intuition lying behind an hybrid functional is 
rather clear. While LDA or GGA calculations strongly 
underestimate the gap, Hartree-Fock calculations overes- 
timates it typically by more than a factor of two. By 
changing the mixing from to 1, one has a continu- 
ous change between local KS and Hartree-Fock, and an 
essentially linear variation between the respective gaps. 
Therefore, to obtain the correct experimental gap one 
simply has to use an appropriate mixing in the functional. 
This value can be determined from a fit to a series of sys- 
tems, and set to around a ~ 0.2 — 0.3. This choice gives 
very good results for a large class of systems, but it usu- 
ally fails when the gap is very large or very small. But 
what is the physical meaning of the mixing parameter? 
To answer this question, it can be instructive to move 
away from a generalized KS picture 14], and consider 
the hybrid as an approximation to the self-energy £ . In 
the GW approximation the latter can be written as 



E(r,r';w) = S sX (r, r') + E r0Bt (r, r'; w) , 



(1) 



with £ s x(t", r') being the statically screened-exchangc 
(sX) term, and E rest (r, r'; u>) containing the static 
Coulomb hole and dynamical contributions. If the 
screening in the sX term is replaced by an effective static 
dielectric constant e,^ = 1/a, and E r cst is modeled by 
the static and local part of the hybrid functional 15[ the 
quasi-particle equation has the same form as t he g eneral- 
ized KS equation solved for hybrid functionals [161 ] . From 
these arguments we can conclude that the physical value 
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FIG. 1. (Color online) Inverse of the dielectric constant cal- 
culated with the PBE as a function of the optimal mixing 
parameter. The line a pt = l/ e ™ E is a guide to the eye. 



for the mixing parameter a is related to the inverse of 
the dielectric constant of the material at hand. Such a 
link has been suggested to explain also variations of band 
gaps with respect to small structural changes [18j ]. 

In Fig. [T] we show that, for a large range of different 
materials, the value of l/e^ - obtained with ABINIT [l9| 
within the PBE approximation - is approximately pro- 
portional to the optimal mixing parameter awt- The 
latter, calculated with the computer code VASP [10(, is de- 
fined as the fraction a of Fock exchange of a PBEO [8| hy- 
brid functional that reproduces the experimental band- 
gap of the material. The correlation is evident, despite 
the fact that the DFT-PBE calculations systematically 
overestimate the dielectric constants with respect to ex- 
periment. 

In Fig. [5] we compare the band-gaps calculated using 
a = 1/ e™ E as the mixing parameter of a modified PBEO 
hybrid functional (red stars labeled PBEOeoo), compared 
to experimental data and other theoretical results. A 
more detailed comparison together with average errors 
can be found in the table given as supplemental mate- 
rials. We can see that the hybrid PBEO (a = 1/4) [1] 
already improves dramatically the result with respect to 
a local PBE functional, bringing the calculated gaps to- 
wards the experimental values with an average error of 
less than a 30%. However, this number alone hides the 
fact that PBEO gives excellent results for some interme- 
diate band gap materials, like diamond, BN, A1N, etc., 
whereas it fails both for small band gap materials (like 
Si, or Ge), overestimating their gaps, and for large band 
gap materials (like the rare-gases) where the gaps are un- 
derestimated. This is not surprising: In materials like Si, 
electrons are delocalized and easily polarizable, leading to 
a strong screening and small mixing (the optimal mixing 
is actually ev op t = 0.12). For Ne, electrons are localized, 



and screening is basically nonexistent (a pt = 0.70). This 
effect is captured by the dependence on the dielectric con- 
stant of the mixing with a = 1/eoo- Such a calculation 
decreases overall the error by almost a factor of two. The 
remaining error is dominated by the large underestima- 
tion of the gap for materials like Si, Ge, or GaAs. This 
point to the fact that, in order to predict good gaps ap- 
proaching the metallic limit, one needs a finite amount 
of Fock exchange that is not accounted for by the sim- 
ple 1/eoo model for the mixing, as it can also be seen in 
Fig.ffl 

The eoo-dependent mixing is hence physically moti- 
vated and can yield good band gaps, however, the pro- 
cedure to obtain it (i.e. the calculation of the dielectric 
constant) is fairly cumbersome and therefore often un- 
practical. The best option would be to find an estimator 
for the dielectric constant from quantities readily avail- 
able from the ground-state. To obtain such relation, we 
first observe that there is a strong correlation between 
the electronic dielectric constant of the material and the 
energy gap pcj . Of course, it is not desirable to have a 
functional depending explicitly on the band-gap of the 
material, so the question is if one can find an estimator 
of this quantity based solely on reduced densities. 

In fact, several density- functional estimators of a "lo- 
cal gap" have been proposed in the past years. For ex- 
ample, Gutle et al. [2l( proposed to use the quantity 
G = ||Vn| 2 /n 2 to define the gap locally. Their argu- 
ments for the use of this quantity were based on the 
asymptotic expansion of the function G for finite sys- 
tems (G will reduce to the ionization energy) and on 
the observed piecewise exponential behavior of the den- 
sity [22l |. More recently, the quantity |Vn|/n was used 
to model a position-dependent scree ning function in a 
so-called "local- hybrid" functional IH 24 1 . Also the von 

|Vn| 2 /8n has 



Weizsacker kinetic energy density tw 
been used to define a "local band-gap", and inserted 
into a "local-hybrid" functional that turns off the exact- 
exchange term when this local gap has metallic charac- 
ter [H. 

To obtain a global estimator of the band-gap of the 
material, and therefore of its static dielectric constant, we 
can average the local estimator over the Brillouin zone. 
We will follow the idea contained in the meta-GGA of 
Tran and Blaha (TB09) [26| and define the quantity 
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where the integral is over the unit cell of volume V^oii- 
We note that the quantity g is very similar to the av- 
erage involved in the calculation of the ctbo9 parameter 
of TB09, and is quite stable regardless on which theory 
is used to evaluate the density. In fact, we verified that 
using as input either a PBE or a Hartree-Fock density 
leads to only minor differences in its value. Our basic 
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FIG. 2. (Color online) Electronic band-gaps calculated as differences of generalized KS eigenvalues for a series of semiconductors. 
All calculations were performed at the experimental lattice constant. The symbols labeled PBE0 Eoo , PBEO m i x and HSE06 m ix 
are the results obtained with the hybrid functionals proposed in this Article. TB09 results come from Ref. [2r3 . 



hypothesis is that the mixing parameter can be written 
as a function of the parameter g. 

To proceed we need to specify the local part of the 
hybrid. We chose to use PBEO form Q which is the basis 
for the screened hybrid Heyd, Scuseria and Ernzcrhof 
(HSE) 13]. As we expected, there is a clear correlation 
between the value of a op t and g; it can be quite well fitted 
by the simple function a(g) 

a = -1.00778 + 1.10507 g, (3) 

with g in atomic units. 

Analyzing the resulting gaps displayed in Fig. [5] (la- 
beled PBE0 m i X ) we realize that fixing the mixing param- 
eter according to Eq. [3] 23] reduces the mean average er- 
ror to slightly more than 14%, much better than PBEO, 
and slightly better than the HSE06 [3l| screened hybrid 
functional. Furthermore, in contrast with PBEO and 
HSE06, our density-dependent mixing describes equally 
well small, medium and large gap systems. The largest 
errors arise for d-electron materials like ZnO where our 
recipe overestimates the mixing parameter and therefore 
the electronic gap. This is due to the fact that local- 
ized d-states give too large contributions to the mixing 
through strong density variations. A possible solution 
would be finding a more pertinent density-estimator for 
those d states. 

In view of the success of screened hybrids in improv- 
ing the accuracy of PBEO, we applied our construction 
also to the HSE06 functional [3l[ . In this case, the phys- 
ical interpretation of the mixing parameter as an inverse 
screening is considerably more complicated, as screening 
is already present to some extent in the range separation. 
Following our protocol, we arrived at the following fit for 
the mixing parameter 

a = 0.121983 + 0.130711 g 4 . (4) 



We remark the different power of g in the expression. 
This is due to the fact that the screening already present 
in HSE06 decreases considerably the strength of the 
Hartree-Fock term, increasing the values of a pt required 
to reproduce the experimental band gap of small-gap sys- 
tems. This is actually the cause of the success of HSE06: 
for example, for Si a op t is now 0.24, which is very close 
to the actual mixing of HSE06 (a = 1/4). Results for 
band-gaps using HSE06 and our new mixing scheme (la- 
beled HSE06 m j x ) are shown in Fig. [2]). Our new mixing 
scheme brings down the HSE06 error from around 17% 
to 10%, of the same order of magnitude as the error in- 
curred by the GqWq approximation and the new TB09 
meta-GGA 0. 

We want now to compare our approach to the the 
new TB09 meta-GGA. The physical interpretation of our 
functional implies that the mixing parameter should take 
values between and 1 , while in TB09 the corresponding 
parameter Ctbo9 is always larger than one. This differ- 
ence stems from the fact that the TB09 functional is a 
purely local potential. Therefore, the band-gap defined 
in terms of total energy differences should be equal to 
the difference of the KS eigenvalues plus the derivative 
discontinuity of the xc potential. It is known that the 
Becke-Johnson potential (2?J (upon which the TB09 func- 
tional has been constructed) reproduces to a very good 
extent the derivative discontinuity of exact exchange for 
molecular chains (28[. Therefore, it is reasonable to ex- 
pect that the DFT band-gap with the TB09 functional, 
after adding the derivative discontinuity, would actually 
become much larger than experiment, and possibly even 
larger than the Hartree-Fock gap since ctbo9 is always 
larger than one, in agreement with what was proved by 
Griming et al. [29j. Of course, the aim of TB09 is to ob- 
tain the gaps simply as eigenvalue differences, and there- 
fore require ctbo9 > 1 in order to compensate for the un- 
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derestimation of the eigenvalue gap in the Becke- Johnson 
functional [30j |. In this sense, it is similar in spirit to the 
pragmatic xa approach (32j . 

From Fig.[3]it is clear that the errors of TB09 is compa- 
rable to ours. On the one hand, the meta-GGA is clearly 
much lighter from the computational point of view than 
any hybrid functional. On the other hand, TB09 is an 
approximation for the xc potential, and is thus incapable 
of yielding total energies (in contrast to our approach) . 
In fact, it can be proved that TB09 is not the functional 
derivative of any energy functional, and therefore violates 
serious constraints, like the zero- force theorem (33|. As 
a consequence, such functionals do not allow to calculate 
structural properties. Our approach, instead, can also 
be used to calculate total energies and structural prop- 
erties. We tested our functionals and found that they 
give relaxed geometries as good as the standard PBEO 
and HSE06 (with lattice constants better than 0.7% for 
the cubic semiconductors considered here). Moreover, we 
mention that the TB09 meta-GGA also inherited some 
of the problems of the Becke- Johnson functional [i3| on 
which it is based (for example, TB09 is not size-consistent 
and gauge invariant [HI). 

It is clear that the averaging procedure in Eq. (0) is 
only meaningful for bulk systems and will fail for finite 
systems, slabs, interfaces, etc. This issue can be fixed 
by converting the global g (and thus the uniform screen- 
ing) into a local function function <?(r), i.e. transforming 
the functional in a local hybrid. The function g(r) can 
be determined by restricting the integral of Eq. ^ to a 
neighborhood of r, with a range related to the screening 
length of the system. This screening length can be taken 
as constant, as in the HSE functionals, or can be a local 
functional of the density (but clearly not a function of 
the average density). This procedure would not change 
significantly the results for the systems studied here as 
they all have very small unit cells. 

In conclusion, we proposed a scheme to calculate on- 
the-fly the mixing parameter in hybrid functionals de- 
pending on the density of the system. In this way, the 
average error on the values of the energy gap were con- 
siderably reduced with respect to the original hybrid 
functionals. The resulting band-gaps are roughly of the 
same quality as those obtained using a GW approach or 
the new meta-GGA of Tran and Blaha. Moreover, this 
method assures also an excellent description of the struc- 
tural properties. These improvements are obtained with 
no increase of computational time with respect to a fixed 
mixing parameter. 
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